#!/usr/bin/perl
# 2012/12/14 Making association map based on the distance curve

use strict;
use warnings;
use IO::File;
use Getopt::Std;
use Carp qw(croak);
$| = 0;

if(@ARGV != 6 or $ARGV[0] eq '--help'){
	die "Usage : $0 -h [resolution conversion table] -o [output association file] -p [distance normalize parameter]\n";
}

my %opt;
getopts("i:h:o:p:", \%opt);
my $FILE_conversion = $opt{h};
my $FILE_out = $opt{o};
my $FILE_para = $opt{p};


#---------------------------------------
# read parameters
#---------------------------------------
my %Param;
my $fh_param = IO::File->new($FILE_para) or die "cannot open $FILE_para: $!";
while($_ = $fh_param->getline()){
	s/\r?\n//;
	my ($type, $key, $value) = split /\t/;
	$Param{$type}{$key} = $value;
}
$fh_param->close();

sub PredictValue{
	my ($type, $x) = @_;
	if($type eq 'inter'){
		return $Param{'inter'}{$x};
	}else{
		my $value = $Param{$type}{'a'} + $Param{$type}{'b'}*exp(-1 * $Param{$type}{'c'}*$x) + $Param{$type}{'d'}*exp(-1 * $Param{$type}{'e'}*$x);
		return $value;
	}
}

#---------------------------------------
# read conversion table
#---------------------------------------
my @binLists;
my $fh_conv = IO::File->new($FILE_conversion) or die "cannot open $FILE_conversion: $!";
while($_ = $fh_conv->getline()){
	s/\r?\n//;
	my ($id, @hin) = split /\t/;
	my ($chr, $middle);
	if($id =~ m/(chr\d):(\d+)\-(\d+)/){
		$chr = $1;
		my $start = $2;
		my $end = $3;
		$middle = int(($start + $end) / 2);
	}else{
		die "bin name is something wrong: $_\n";
	}
	push @binLists, [$chr, $middle];
}
$fh_conv->close();


#---------------------------------------
# output
#---------------------------------------
my $fh_out = IO::File->new($FILE_out, 'w') or die "cannot write $FILE_out: $!";
for(my $i = 0; $i < @binLists; $i++){
	my ($chr1, $pos1) = @{$binLists[$i]};
	for(my $j = 0; $j < @binLists; $j++){
		my ($chr2, $pos2) = @{$binLists[$j]};

		my $value = 0;
		if($chr1 ne $chr2){
			my $type = '';
			if($chr1 eq 'chr1'){
				if($chr2 eq 'chr2'){
					$type = 'chr1_chr2';
				}elsif($chr2 eq 'chr3'){
					$type = 'chr1_chr3';
				}
			}elsif($chr1 eq 'chr2'){
				if($chr2 eq 'chr1'){
					$type = 'chr1_chr2';
				}elsif($chr2 eq 'chr3'){
					$type = 'chr2_chr3';
				}
			}elsif($chr1 eq 'chr3'){
				if($chr2 eq 'chr1'){
					$type = 'chr1_chr3';
				}elsif($chr2 eq 'chr2'){
					$type = 'chr2_chr3';
				}
			}
			$value = &PredictValue('inter', $type);
		}else{
			my $distance = abs($pos1 - $pos2);
			$value = &PredictValue($chr1, $distance);
		}
		$fh_out->print("$value\n");
	}
}
$fh_out->close();




